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It is demonstrated that the commonly applied self interaction correction (SIC) used in density 
' functional theory does not remove all self interaction. We present as an alternative a novel method 

' which, by construction, is totally free from self interaction. The method has the correct asymptotic 

l/r dependence. We apply the new theory to localized /-electrons in praseodymium and compare 
5-H ' with the old version of SIC, the local density approximation (LDA) and with an atomic Hartree-Fock 

calculation. The results show a lowering of the / level, a contraction of the / electron cloud and a 
lowering of the total energy by 13 eV per 4/electron compared to LDA. The equilibrium volume of 
the new SIC method is close to the ones given by LDA and the older SIC method, and is in good 
agreement with experiment. The experimental cohesive energy is in better agreement using the new 
SIC method, both compared to LDA and another SIC method. 
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cy:j ; INTRODUCTION 
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^ Density functional theory (DET) states that the ground state energy of a many particle system is an exact functional 
^ ; ' of the total particle density pa. Unfortunately the exact functional is unknown and in order to have a practical mean 
, of calculating the ground stat&.energy one maps the interacting particle system onto a non-interacting, for which the 
^ ■ total energy can be written aa3 

: E[p\^T[p] + [p] + Exc [p] + ^NN + ^^cN [p] , (1) 

« ; 

O ' where T[p] is the one-electron kinetic energy (given by ^niei — JVcsp)i £^h[p] is the Hartree interaction. Exc[p]i 
I the exchange-correlation (XC) potential, is supposed to correct for the approximations introduced when replacing 

the true many-body kinetic energy by T[p\ and the correct electron-electron interaction by -BHip]- In addition -Enn 
is the electrostatic nucleus-nucleus interaction, and _Bcn[p] is the electron-nucleus interaction. Since the exact XC 
^ potential is not known, one has to rely to approximations. The simplest (and most frequently used) approximation 
, is the LDA, which is a parameterized XC functional using a homogeneous electron densitjf- Approximate expressions 
' for the total energy combined with efficient methods for solving the Kohn-Sham equationscl have resulted in powerful 
, methods for studying many materials properties such as magnetism, chemical bonding, hardness, elasticity as well as 
surface related phenomena, including relaxation, reconstruction and chemisorption. Quite frequently LDA calculations 
T-H reproduce experimental results with an error smaller than a few percent. However, notable exceptions exist. One 
problem is that the approximations leading to the LDA introduce an unphysical self interaction (SI) of a particle with 
itselL For instance, in the local density approximation the Hartree part contains the mean-field interaction energy of 
an electron with itself. An exact XC functional would have a SI part which would exactly cancel the SI of the Hartree 
^ ' term. The result would be a theory totally free from self interaction, giving the exact ground state energy. It may be 
r I noted here that for systems with broad bands, where the electron states form Bloch states, self-interaction becomes 
• neglibly small, whereas for systems where the electrons are localized or close to localization the self-interaction problem 
Q , becomes serious. A self interaction correction (SIC) method ofi-jBFT was presented in 1981 by Perdew and Zungercl 
O ■ (here referred to as PZ). Although there exist improved schemeaaS the PZ-method is the method predominantly used. 
One of the virtues of the correction of PZ is that for the exact XC functional the added correction becomes zer^q.r^he 
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method of PZ has been applied with some success to difEer|eat systems, e.g. localized electrons in Ce and PvB lla 



1-D Hubbard modelt3 and transition metal mono oxidesQ^cZHiS. Although the PZ version of correcting for SI has been 
^ . successful in many applications we question its ability to remove all self interaction and suggest, as an alternative, a 
- - ' novel way of performing LDA-SIC calculations, which is totally free from self interaction. 

EXCHANGE AND CORRELATION IN SELF INTERACTION CORRECTIONS 

The self interaction corrected XC energy, within the approach by Perdew and ZungerO, applied to LDA, is normally 
written as 
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E'^§{n]=EMF[n]-Y: 



(2) 



The basic idea behind this equation is to subtract the interaction energy of an electron with itself from the total 
energy. The second term on the right hand side of Eqn. (j^) corrects for the self-Hartree interaction but we choose to 
include this correction to the XC functional since formally one may introduce all terms there which are not found in 
the kinetic energy functional and the Hartree functional. The last term on the right hand side of Eqn.(||) is called the 
self XC. The in exc(pf (r),0) indicates that the correction is made separately for each spin-channel. Minimization of 
the SI corrected energy functional with respect to the density leads to Kohn-Sham like equations with a SI corrected 
effective potential, Vcff- The corrections coming from SI are normally referred to as the SIC part of the potential. The 
correction for the self Hartree interaction is exact, but as will be argued below the correction for self XC is imperfect. 
To illuminate the problem, we consider an LDA functional (Eqn.(Q)) and take a closer look at the XC functional. 
It is written as i?xc(p(r),0) = / d'^rexc(p(r), O)p(r), where one term in the XC energy density, exc, is proportional 
to p^/"^. Let us now write the total density as ptot — P + pf , where pi is the density given by the orbital which we 
are interested to correct for self interaction, p is the summed density of all other orbitals and ptot the total electron 
density. It is now easy to see that the correction of PZ, as defined in Eqn.(H), does not remove all self XC, since with 
the above expressions for ptot,p and pf , we obtain from Eqn.(||) an expression for the exchange part of the energy 
functional 



dV(p+pn^/^ 



iPtr^'- 



(3) 



Eqn.(||) is not equal to j cPr{p)'^/^ and hence does not represent a contribution to the energy functional that is free 
from self interaction. 

In a similar way the correction of PZ does not remove the self XC interaction in the effective potential of the 
Kptn-Sham equations. The XC part of the effective potential, Vxc^ of the Kohn-Sham equations is in LDA written 
all 



(4) 



where it may be seen that the exchange potential (^[i^{vs)1(ij^ contains a term proportional to p(r)^/'^ (for a description 
of terms in Eqn(^,see Ref.Ej). The PZ corrected effective potential for electron state i involves, analogous to Eqn.(|), 
the following terms, 



PZ 



p(r') 



-St' + A(r)([p(r) + pf (r)]!/^ - pf (r) ^3) 



(5) 



and we observe that due to the non-linear dependence on the electron density, y^^^ is dependent on pi and hence SI 
effects are present in the effective potentials well. This was also pointed out by other authors 



A NOVEL METHOD OF SELF INTERACTION 



We have now pointed out why the SIC form of PZ does not remove all the SI and below we suggest a method which 
does this. In general, the approach is to identify all contributions in the LDA energy expression (Eqn.(|l|)) and the 
effective potential which contain SI and to remove these interactions. Writing ptot = P + pf : one may obtain an energy 
expression which is free from self-interaction, by applying the correction to the densities instead of the energies. The 
energies for pf are calculated from p instead of the total density ptot- LDA contains a term £'[p-f pf](p-|-pf) whereas 
the wanted expression is E\p^ p'[]p + E{p)pf [E denoting the total energy functional) if we only apply SIC to orbital 
pf . The expression for the total energy, when formulated as — E 



LDA 



E'' 



IS 



SIC 



^E, 



LDA 



SICn 
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Pf(r)pf(r') 



r r 



[vM-vM]p1d\, 



(6) 



where E^^j^ 



gSiC^ is the normal LDA total energy expression (Eqn.(|i|)) evaluated with a KS-eigenvalue obtained from 
Eqn.(|^) below, the effective potential. To arrive at the expression in Eqn.(^) we have explicitly subtracted the parts 
in the LDA expression which contain SI and thereafter added a SI free part. 
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The effective potential to be used in the KS-equation is constructed in the usual way 



V^ff (p, ft , r) = + ^'^N, (7) 

where V^n is the electron- nucleus potential 2Z/r. Explicitly taking the functional derivative, we get the potential 

P^^"> ^3 / _ A/ X-/ xl/3 



I 

This effective potential is by construction totally free from SI if one uses the electron density p in the functional 
A(r). The Hartree part in this equation coincides exactly with the Hartree part in Eqn.(^). In short, the potential 
for the orbital pf can be written as V{ptot — Pf,^) = V{p,^), instead of the normal orbital independent potential 
V{ptot,^), or the PZ-SIC: V{ptot, r) — V{pf , r). This procedure is applicable to all approximations for the XC energy 
but here, for simplicity, we choose to discuss the LDA. In the asymptotic limit r ^ oo the effective potential has the 
correct - dependence, since we subtract one electron from the screening cloud of the neutral atom. One can arrive 
at Eqn.(^) by starting from the Hartree- Fock expression for the effective potential, which is free from self interaction, 
and then replace the Fock part by A(r)yo(r)-'^/'^. This replacement then hopefully both simplifies the effective potential 
by making it local as well as includes correlation effects. 



NUMERICAL TESTS ON THE METHOD 



We have implemented this new method, and the method by RZ., to be able to make a comparison between the two, 
in a full potential linear muffin-tin orbital (F|E,-LMTO) methadpL Since Pr is a material with 2 localized f-electrons 
which has been studied extensively by LDAOOtJ, PZ-SICll3il3 and the orbital polarization methodlj we choose 
to test our theory for this element. For simplicity Pr was put in an fee structure (observed structure is the dhcp 
structure), and we used 2197 k-points in the whole Brillouin zone. The SI corrected functionals were implemented 
only for the 4f core electrons. In the implementation we have ignored thejion-orthogonality of the orbitals introduced 
by SIC since in atomic like problems it has been shown to be negligibleo. The results for praseodymium are shown 
in the Table where the total energy contributions are listed, and in Figures ^ and ^ where the electron density and 
different contributions to the effective potential are displayed. 



Quantity 


LDA 


PZ-SIC 


SIC 


4/ KS-eigenvalue'' (Ry) 


-0.054167 


-0.771376 


-1.537387 


Equilibrium volume'^ {Rsws in a.u.) 


3.826 


3.842 


3.916 


Correction to energies 








Kinetic (Ry) 




-0.1239143 


-1.0536124 


Hartree (Ry) 




-1.8285461 


-2.0564759 


XC (Ry) 




1.7158233 


0.522432025 


Total Energy (Ry) 


-18'471.575 781 


-18'471.825 455 


-18'473.480 380 



TABLE I. Comparison between LDA, PZ-SIC and the current SIC method. The energies come from the corresponding 
minima in volume. The energies are the corrections to the (LDA) 4/-electron energies given by Eqn.(^ for two 4f electrons. 
''SIC From Ref.i. 

''With respect to the Fermi energy. 

'^The experimental equilibrium volume is 3.818 a.u. 



3 



2.5 



2.0 



1.5 



13 



g 10 



0.5 



0.0 



LDA 
PZ SIC 
SIC 

HF (atom) 




0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 

r (a.u.) 

FIG. 1. The 4/ radial density using different approximations for exchange and correlation. 
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FIG. 2. Fig. A shows the difFerent contributions to the SIC potential in the PZ-method. Fig. B shows the corresponding 
potentials for the current implementation of SIC. 

By inspecting the 4/-density (shown in Figjl]) one observes that the orbitals become more contracted when SIC is 
used. Notable is also that the radial extent of the orbital in the new SIC method is largely reduced compared to the 
PZ-SIC method. We also note that compared to an atomic Hartree-Fock density the /-cloud is actually even more 
contracted. From Fig.Hb one sees that for the current SIC method the main contribution to the potential entering 
the KS equation comes from the Hartree part, the XC part provides a very small contribution. This can be compared 
with the PZ-SIC potential shown in Fig.^ja, where the same Hartree part is used, but there is a big correction to 
the XC potential. The reason why the correction to the XC potential in this implementation is so small is that the 
potential has the form Vxc — ^c(Ptot — pf) instead of the LDA Vxc = ^c(Ptot) and for large values of ptot the 
difference between them, ptot — Pi = P ~ Ptot, is very small. As a side remark we note that the electronic structure 
and density of states of the valence states is almost unchanged (no SIC is applied to the valent electrons), showing 
that the contraction of the 4/cloud due to SIC only marginally influences the energy-band dispersion of the valence 
states. 

From the Table it is seen that the KS eigenvalue is significantly lower in the new SIC method compared to LDA 
and PZ-SIC. The fact that the /-orbital localizes more gives rise to a larger screening for the valence electrons, which 
as a result become more extended and for this reason a larger equilibrium volume is found. However, as is seen from 
the table the difference in equilibrium volumes between the different approximations is only a few percent and they 
agree rather well with experiment. The small deviation from expepment can be explained by a contribution to the 
cohesive energy from delocalizedJ'-electrons as described previouslyEj. The cohesive energy was previously calculated 
from LDA to be 5.09 eV/atomE3, a value which is 0.55 eV/atom higher than experiment. Our calculation lowers 
the cohesive energy with 0.29 eV/atom compared to the LDA calculation, and is hence in better agreement with 
experiment. The PZ calculation lowers the cohesive energy by only with 0.17 eV/atom compared to LDA. 

In the Table we also note that the self-Hartree part is larger in our method of SIC than in the PZ-method, this 
is due to the more contracted 4/-orbital. In the PZ-SIC the total energy changes only by 0.25 Ry since the Hartree 
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correction is balanced by a change in the XC energy. In our implementation of SIC the Hartree energy is not balanced 
to the same extent by the XC energy, thus giving a change in total energy by 1.90 Ry. We also compared the 4/- 
eigenvalue with an atomic, relativistic Hartree-Fock calculation on Pr'^^ with a 4? configuration which is free from 
SI (but without correlation effects). This calculation gives a value of e^f =-2.870 Ry, which may be compared to an 
atomic-like calculation (a band-calculation with lattice expanded) which gives -2.3062 Ry for the new SIC method, 
-1.175 Ry for the PZ-SIC method and -0.444 Ry for a normal LDA calculation. Hence the current SIC method is 
closer to a Hartree Fock description than LDA or PZ-SIC. This finding is consistent with the discussion above that 
Eqn.(|^) may be viewed as an approximation of the Hartree-Fock equation. 

The present SIC method lowers the total energy compared to the method of Perdew and Zunger. For this reason it 
may not be as suitable for describing localization-delocaiisation transitions in the rare-earth metals, since the PZ-SIC 
method reproduces the Mott transition in Ce very welE'EI. However, since the PZ-SIC method does not remove all 
self interaction it would seem that the good description of this method in describing the localization-delocalization 
transitions in the rare-earths is somewhat fortuitous. 

We also remark that as long as a Hartree term (which depends linearly on the electron density) is introduced in 
the energy functional, any local XC functional exc(/3i(r), 0) which is a functional of the total density must also have 
a linear dependence on the electron density, if SI is to be absent in the functional. 

By constructipn the LDA XC potential cancels some of the unphysical self interaction, e.g. in hydrogen 95% of the 
self interactionEI is cancelled. Application of the presented SIC method to the hydrogen-atom would yield an exact 
result (as would the PZ-SIC). Since the application of PZ-SIC normally improves the LDA result, we believe that the 
currently described method of SIC will be successful when applied to systems like the 1-D Hubbard model and the 
bandgap problem in NiO, since it properly corrects for all self interaction. It may also be a good starting point for 
finding an exact XC functional of localized or nearly localized electron systems. 



CONCLUSION 



In conclusion, we have presented a novel way to perform SIC-corrections in DFT calculations, and compared it to an 
older, not fully SI free, method. In its form the suggested correction is similar ta the corrections of the non-spherical 
4f density proposed by Ilrooks et al. in calculations of crystal field excitatiDnsE3. The method presented is closer to 
the suggestion of NesbetQ and differs from the method of e. g. Manoli et acA whose approximation is more similar to 
that of Ref.Q. We have applied the new SIC method to praseodymium metal and showed that it results in stronger 
localization of the f-shell, thus correcting more for self interaction. Despite their large differences in describing the 4f 
electrons, the calculated equilibrium volumes of LDA, PZ-SIC and the presently proposed SIC method are all quite 
similar and agree with the experimental values. This is a result of that the valence states (that are not calculated 
using any SIC potential) are rather independent of the shape of the 4/-density and it is namely the valence states 
that are responsible for the chemical bonding. On the other hand, the cohesive energy is better described by the new 
SIC method, both compared to LDA and to PZ-SIC. 
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